clear
clc
tic

% parameters
np = 5000; % increases the pixel-density of the graph
p1 = linspace(0,1,np); % recognition probability of player 1
p2 = linspace(0,1,np); % recognition probability of player 1

sim = nan(np,np);  % points on a simplex

%conditions 

for i=1:np
    for j=1:np
        if p1(i) + p2(j) <=1 % probabilities must add up to 1
            sim(i,j)=1;
            if p2(j) > (p1(i)^2)/(1-p1(i)) % V2(p) > V1(p)
                sim(i,j) = 1;
                if p2(j)<(1-(p1(i)^2))/(2+p1(i)) % V3(p) > V2(p)
                    sim(i,j)=4; % This is equal to 4 only when V3(p) > V2(p) > V1(p)
                end
            elseif p2(j)<(1-(p1(i)^2))/(2+p1(i)) % V3(p) > V2(p)
                sim(i,j)=1;
            end
        end

    end
end

plot1 = nan(1,np);
plot2 = nan(1,np);
plot3 = nan(1,np);
plot4 = nan(1,np);
plot5 = nan(1,np);
plot6 = nan(1,np);
plot2a = nan(1,np);
plot3a = nan(1,np);
plot4a = zeros(1,np);
plot7 = nan(1,np);

for i=1:np
    plot1(i) = 1-p1(i); % boundary of the simplex
    plot2(i) = max(0,min(1-p1(i),(p1(i)^2)/(1-p1(i)))); % above this line, V2>V1
    plot3(i) = max(0,min(1-p1(i),(1-(p1(i)^2))/(2+p1(i)))); % above this line, V3>v2
    plot4(i) = min(plot2(i),plot3(i));
    plot5(i) = max(0,min(1-p1(i),p1(i))); % below this line, p1>p2
    plot7(i) = max(0,min(1,(0.5*(1-p1(i))))); % above this line, p2>p3
    if plot3(i)>=plot2(i)
        plot2a(i) = plot2(i);
        plot3a(i) = plot3(i);
    end
    if p1(i)<=0.434
        if p1(i)>=0.382
            plot4a(i) = plot4(i);
        end
    end
end





plot(p1,plot1,'black')
hold on
patch([1/3,.366,.434,.414,1/3],[1/3,.366,1/3,.293,1/3],'black','FaceAlpha',.4)
xlabel('Recognition probability of player 1: p_1')
ylabel('Recognition probability of player 2: p_2')

%% Explanation of points of intersection
% (1/3,1/3): p1 = p2, p2 = p3
% (.366, .366): p1 = p2, V2(p) = V3(p)
% (.434, 1/3): V1(p) = V2(p), V2(p) = V3(p)
% (.414, .293): p2 = p3, V1(p) = V2(p)
% 

%%
% uncomment below to see different plots, you might want to reduce 'np' to
% 250 though

%plot(p1,plot1, p1, plot2, p1, plot3, p1, plot5, p1, plot7)



toc
